Stripe Formation in Bose-Einstein Condensates with Large Numbers of Vortices 
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We analyze recent experiments at JILA which observe periodic formation of stripes in a rotating 
Bose condensate [Engels et al. Phys. Rev. Lett. 89, 100403 (2002)], and explain this behavior 
as a consequence of kinetic energy driven evolution, where the vortices follows a quadrupolar flow. 
We demonstrate that the structure of vortex lattices in the lowest Landau level is a useful tool for 
understanding such complex phenomena, even when the system is not in the quantum Hall regime. 
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In the last year and a half there have been remark- 
able developments in producing Bose-Einstein conden- 
sates (BEC) with large numbers of vortices Q]. These 
experiments provide a unique setting for studying the 
dynamics of vortex lattices. Such experiments would be 
nearly impossible in a conventional ultra-low tempera- 
ture laboratory. It is conceivable that even higher an- 
gular momentum states can be achieved, eventually pro- 
ducing quantum Hall physics. 

The novelty of the dynamics at high angular momen- 
tum is demonstrated in a recent experiment at JILA^). 
In this work, after a vortex lattice (with rotation fre- 
quency f2) is created, the stirring force is turned off and 
a small quadrupolar deformation is suddenly introduced 
in the trap. The vortex lattice is then rotating against a 
static anisotropic potential. Even for small distortions, 
the atom cloud undergoes significant ellipsoid deforma- 
tion. Surprisingly, the vortex lattice transforms peri- 
odically into a set of parallel stripes, as if the vortices 
"melt" into sheets and reappear in a coherent fashion. 
The stripes appear six times over the rotating period 
T = 2ir/fl, whenever the basis vector of the vortex lattice 
coincide with the minor axis of the cloud. 

There is little question that this intriguing behavior is a 
general feature of rotating Bose gases with large angular 
momentum (L) , and hence a consequence of the unique 
physics in this limit. From the view point of energetics, 
the most significant change as L increases is the rising 
importance of kinetic energy, which ultimately forces the 
system into the lowest Landau level. We quantify this 
consideration by showing that at sufficiently high L, the 
interactions only perturbatively enter the dynamics, and 
the time evolution of the vortex lattice is driven by ki- 
netic energy. We calculate analytically the time evolu- 
tion of the condensate, choosing a variational form of the 
initial vortex lattice based on its equilibrium structure 
at even higher angular momentum (i.e. the mean field 
quantum Hall regime) Q. Our calculation explains the 
origin of all the features observed in ref . J|] . The agree- 
ments between our results and experimental observations 
shows that even if the system is not in the quantum Hall 
regime, the physics of lowest Landau level, which plays 
an important role in our calculation, is useful for un- 
derstanding complex phenomena and provides powerful 
approximation schemes. 



I. Kinetic evolution: Justification and evolution 
of density profile: We are interested the solution of the 
Gross-Pitaevskii (GP) equation 

iW(x,i) = (T + V(x) + .g|v]/(x,i)| 2 )vl/(x,i) (1) 

with the boundary condition ^(x, t = 0) = &(r), where 
T = -(?i 2 /2M)V 2 , $(r) is the equilibrium vortex state of 
a cylindrically potential, ^(r, t) is the condensate wave- 
function after a quadrupolar distortion is added to the 
potential at t — 0, and V(x) is the deformed trap at time 
t > 0. We shall treat the problem as two dimensional, 
since the effect of the third dimension can be obtained 
by Thomas-Fermi approximation as in ref.|]. In this 



case, we have V = \M(uj 2 x 2 + u> 2 y 2 ), 



uj(1 



uj y = uj(1 — e), where e is the quadrupolar distortion. 

Since $(r) is an equilibrium rotating state, it mini- 
mizes the GP functional T — (ho + U—flL z — pN), where 
h = T + V, (U) = (g/2) j]$| 4 , p, is the chemical poten- 
tial and L z = —iM • x x V. Writing $(r) = e lS ( r ) y/p{r), 
where S is the phase function and p is the density, we 
have T = J{|^[VVP] 2 + f (v s - ftz x r) 2 p+ [\M{u 2 - 
fl 2 )r 2 — p]p + \gp 2 } where v s = (h/M)\7S is the super- 
fluid velocity. For equilibrium states with large number 
of vortices, the system has a spatial average v s w flz x r, 
and a density is 



Po~{p~ ^M(uj 2 



n 2 )r 2 )/g. 



(2) 



It is important to note that states with large number 
of vortices can only be achieved when — > u>, which 
means that \M{u? - fl 2 )r 2 « V(r) = \Moj 2 r 2 . As 
a result, the g\^\ 2 term in Eq. (|l|) can be replaced by 
a constant (g| v I / | 2 ~ p) at early times where \P w 
We show below, however, that this replacement remains 
valid over the observed period of stripe formation, and 
the dynamics are governed solely by the kinetic energy, 
leading us to describe the dynamics as kinetic evolution. 

Under kinetic evolution, the shape of the condensate 
is described by the matrix — (xi(t)xj(t)} /[(r 2 ) /2], 
where (sci,^) = (x,y), Xi(t) = e lh ° t x l e~ th ° t , and (..) 
means average with respect to the initial state $. Since 

X(t) = C X X + S X P X /{MUJ X ), y(t) = Cyl/ + SyPy/(M0Jy), 

and p^> w MQz x r^P, where c x — cosa-^i, s x = sina^i, 
we have (xaj) = ^(r 2 )8 ij: (p.pj) = \M 2 VL 2 (r 2 )8 i:j - 
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\M 2 uj 2 {r 2 )5ij, (xp y ) w Affi(i 2 ), and hence 



V = 



c 2 +t 2 

C x ty 



Cyt x 



C X ty Cytx 

c 2 + t 2 
S + L y 
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(3) 



where = SiUj/idi, i = x,y. (Eq. (g) can also be derived 
from the trial wavefunction we use later) . Noting that ini- 
tial density in Eq. (g) is of the form p a = p(l — r 2 /(r 2 )), 
where gp = \M{uj 2 - £l 2 )(r 2 ), Eq. (g) implies that with 
similar approximation, the deformed density is p(r) = 
; p ( 1 — r • V^ 1 ■ r/(r 2 )), where the inverse determi- 
nant comes from number conservation. The deformed 
density can be rewritten as 



P(r) = 



M(co 2 -n 2 ) 



y/X+X- 2 5v /A+A_ V A + 



where A± and n± are the eigenvalues and eigenvectors 
of T>, (A + > A_, n + • n = 0), and r± = r ■ h±. To be 
able to replace gp in Eq. (|l|) by a constant, one needs 
gp « \Muj 2 (r\ + r 2 ), or 



1 - Q/U) « 



(5) 



In Fig. la, we plot the A-product in Eq. (||) and the 
aspect ratio ^/A_/A+ of the ellipsoidal density profile 
Eq. (4) as a function of u>t for e = 0.036. It is clear that 
Eq. (5) remains valid up to 2.5 to 3 times of the original 
period 2tt/uj for rotating states with 1 — VL/uj ~ 10~ 2 . 
One can see from fig. la that even for e = 0.036, the 
aspect ratio reaches 50% when tot ~ 37r[||. In fig. lb, 
we have plotted the orientation of the minor axis n = 
(cosx, sinx), where \ can be found from Eq. (||) to satisfy 
tan2% = 7/(0: — 0). We see that n_ starts at x, and then 
rotate by 45° in half a period and stays around that value 
for some time. Eventually it will rotate back to x ~ 0, 
but only at times much longer than that for which kinetic 
evolution is valid. 

For the rest of the paper we shall stay within the time 
interval where kinetic evolution is valid. Eq. (Q) then 
gives 



f(x,i)= J G(x,x';i)*(x') dx', 



(0) 



where G(x,x';i) = (x\e~ lhot / h \x'), h = p 2 /2M + 
M(u 2 x 2 +w 2 j/ 2 )/2, which can be shown to be 

G(r, r';t) = ^ e (i/2)(r-Q-r+r'.Q.r'-2r.p.r') ) ( 7 ) 



and A = {2T:a 2 Y 1 ^^f-, where a 2 = H/(Muj). 

II. The evolving condensate: For later discussions, 
we represent the basis vector bi, b2 of a vortex lattice as 
complex numbers bi = b^ ■ (x + iy), i — 1,2, and define 
t = 62/^1 = u + iv, where u and v are the real and 
imaginary part of r. The entire lattice will be denoted 



as {6}, where b = l\b\ + £2^2, and £±,£2 are integers; 



and the size of unit cell is 



The choice of basis 



is not unique. A natural choice is the "fundamental set" 
(61,62) containing the two shortest non-parallel vectors 
in the lattice. As well known from Conformal Analysis, 
this set is confined in the "fundamental" region (|r| > 
1, \u\ < 1/2, v > 0), and that for any basis (61, 62) one can 
find an image (61, 62) in the fundamental region by basis 
transformation. Because the way the basis vectors are 
chosen, the fundamental region is best for describing the 
stretching of lattices. To distinguish the initial hexagonal 
lattice from others at later time, we label its quantities 
with a superscript Thus, r° = e l7T / 3 , |6J| = |6f|. 

To allow for an analytic calculation of the condensate 
dynamics we use a variational form of the initial state, 



<|(r) = C\\{z-b°)e 



r 2 /(2d 2 ) 



(9) 



where z — x + iy, C is a normalization constant, The 
lengths d and |6°| are variational parameters obtained by 
minimizing the energy ho + U — QL Z — pN. Eq. ([)]) is 
the form that will emerge in the mean field quantum Hall 
regime ||, (in which case d is no longer be a variational 
parameter but is equal to the trap length a exactly) . The 
minimization process is similar to the one described in 
ref. H . Since our results do not depend on the numerical 
value of the optimal parameters but only their general 
features, we shall only list the relevant features and not 
dwell on the details of the minimization process. These 
features include: 

(a) A condensate with large number of vortices has ro- 
tational frequency Q — *■ to, and d — > a. 

(b) The vortices form a triangular lattice with averaged 
vorticity (2ttH/M)/v° = 2Q, « 2u>, which gives v° = ira 2 , 
and \b1\ 2 = I62I 2 = 2vra 2 /V3. 

(c) The overall density of $(r) Gaussian with a width a 



(° 2 



1 1,) such that a = d — ttv c . Since d 



a. 



we have a S> a. 

(d) $ satisfies the lowest Landau level constraint [(p x + 
ip y ) + Tiz/d 2 ]^ = 0. It is straightforward to show from 
this condition that V is given by Eq. (g) when a 2 > d 2 . 

(e) Apart from a proportionality constant, $(r) is the 
Weierstrass a- function^. 

To proceed with the integration in Eq. (||), we use the 
relation between the Weierstrass er-function and the 9- 
function to write $ as ||. |(J 

$(r) = 0(q o • r, To)e <^y? H 2 <) l{2d>) (1Q) 

where q Q — (x + iy)/b°, t = e 27 ™/ 3 , and 
00 

6((,T)=- Y (-l)"e MT (n+l/2) 2 e 2«C(«+l/2). (n) 



n— — 00 



We have set the normalization constant in Eq. ( ]To| ) to 
1 because it is not important for our discussions. To 



3 



proceed further, we rewrite $ as 
*(r) = 6(q -r,T )e-^ M °- x 7 



Mo = 



-73- ^ -l-^ 



1/0 -i 



2 \ -i 2 



(12) 
(13) 



The approximate form is due to properties (a) and (c) 
mentioned above. 

Because of the Gaussian structure in Eqs. (|IT|), (12), 
and (^) , the integral (|J) is also a Gaussian. After sorting 
out terms containing (n + 1/2) 2 and (n + 1/2) in the 
exponent, ^(r, t) can be recast into a form similar to 
Eq. ©, 



#(r, f) = 66> (q • r, t) e - 
= PA- 1 ^ - iQ, A = X + iQ, 
q = -iPA"^, q = (x + iy)/6f, 
t = r + r), r) = i2nq ■ A -1 ■ q 
B = 27iv4/ A /Det|A|. 



(14) 
(15) 
(16) 
(17) 
(18) 



Equations (|14j) to (|18|) are the key results of this paper. 
An important feature of Eq. ( |l4| ) is that the zeros (i.e. the 
vortices) of \&(r) continue to organize in a lattice during 
the evolution, with basis vectors (bi and b 2 ) given by 



q • bi = 1, q • b 2 = t, 



(19) 



as a consequence of the quasiperiodicity properties, 9(( + 
l,r) = -0(C,t), and 6{( + r, r) = -9{(, T ) e -^ T + 2 ^ . 
One can see from Eq. ( |l4|) that the vortices fol- 
low a velocity field V(x) = QijXj, where Qij(t) = 
[(9 f bii)Kij + (9 f b i2 )K 2 j] /(27r), and K^s are the recip- 



rocal basis vectors, K, • b , 



2nd, 



■j X j. The symmetric and 
antisymmetric part of Qij correspond to a quadrupolar 
and a rotational flow. Thus, apart from an overall rota- 
tion, the vortices follow the streamlines of a quadrupolar 
flow. 



By properties (a) and (b), the parameters of Eq. (14), 
are approximately d s» a,J6°| 2 w 2na 2 /\/3, which fully 
specifies q, r, M. and V m. While our results hold for 
arbitrary e, we shall focus on small deformations e << 1 
in order to compare with current experiments Q. In Fig. 
2, we have plotted the density of the system at 8 equal 
time intervals between cut = and 6 for e = 0.036 and 
fl = 0.92cj. Stripes appear periodically, becoming more 
distinct at each subsequent appearance. 

To illustrate the relation between lattice deformation 
and stripe formation, wc find the images (bi, 6 2 ) of basis 
vectors obtained from Eq. (^) in the fundamental region. 
The stretching of the lattice is then describe by the ratio 
t = b± I hi =u + iv. Large stretching is given by large | r| 
or by large v since \u\ < 1/2. In Fig. 3, we have plotted 
v and the angle (j> between b\ and x axis as a function 
of tot at the lower panel. We find that (i) stripes are 
most pronounced when v is a maximum, which occurs 
six times within a period cut. Moreover, the sharpness 
of stripe density increases with time; (ii) the maximum 



occurs when <f> is 45°, which is the direction of the minor 
axis after half of the period as shown in Fig. lb. These 
two detailed features are in exact agreement with those 
observed in ref . Q . The discontinuities in (f> are due to the 
mapping of the basis vectors back into the fundamental 
region when the lattice is stretched. The orientation of 
the second basis vector 6 2 is given in the top panel of 
fig. 3, where we plot the relative angle k between b\ and 
b 2 - Once again, the discontinuities of this angle is due 
to the mapping back to the fundamental region. To il- 
lustrate the continuous evolution of the lattice and the 
discontinuous change of the fundamental basis, we have 
plotted the lattices and indicated their fundamental basis 
at different times in Fig. 4. It is also interesting to note 
that at the cusps in fig. 3, the unit cell of vortex lattice 
is close to a square. 

(Ill) The origin of stripe formation: Figure 3 and 
4 show that vortex stripes can be produced by rather 
moderate changes in the basis vectors. This sensitivity 
is due to the special coherent properties of the vortex 
lattice, as revealed by studying the periodic part of the 
density which is calculated from B, 



1 E 



5Ke 



K r 



K 



(20) 



where K =miKi+m 2 K 2 are the reciprocal lattice vec- 
tors, mi are integers, and <?k = (— l) m i+ TO 2+ m i m2 
e -«c|K| / 87r ^/v c /2. Because of the rapid decay of gn with 
increasing |K|, only a few K vectors around the origin 
have non-trivial contribution to the density. For a trian- 
gular lattice, it is natural to organize the K — points into 
concentric hexagonal rings, corresponding to reciprocal 
lattice vectors of roughly equal length, and to number 
these rings, starting from the origin. The same label- 
ing works for distorted lattices, and for sufficiently small 
distortions the first ring, made up of the six shortest 
wavevectors, has the largest contribution to the density. 
As |t| = |6 2 /6i| increases, b\ becomes the shortest vector 
in real space, which in turn makes ±Ki the shortest vec- 
tor in k-space, and this wave-vector dominates the con- 
tribution from the first ring. Similarly, ±2Ki dominates 
the contribution from the second ring, and so on. As a 
result, the density is increasingly dominated by the sum 
of 3k 's containing only integer multiples of Ki, which 
gives rise to a sequence of stripes in real space. 

Final Remarks: We have shown that the mysterious 
periodic stripe formation in ref. Q reflects the kinetic evo- 
lution of the vortex lattice, where vortices remained orga- 
nized in lattice, with each vortex following the streamline 
of a quadrupolar flow. The six occurrences of stripes dur- 
ing a period is caused by an almost periodic deformation 
of the vortex lattice and not the melting of vortices into 
sheets. The method presented here should also be appli- 
cable to other to dynamical phenomena at high angular 
momenta and to more complex condensates [pj . 
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FIG. 2: Density profile of the evolving condensate from Lot — 
to 6 at equal time interval. Dark color means low density. 
Panel (h) corresponds to point (c) on Fig. H. 
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FIG. 1: (a) The product X_ 2 \V 2 and the aspect ratio ^/\- 
for e = 0.036. (b) The angle x between the minor axis n_ 
and the a;-axis. 



FIG. 3: The lower panel shows v (dotted line) and the angle 
4> between b\ and the x-axis as a function of tot. The peaks of 
v correspond to the appearance of stripes, and coincide with 
4> = 7r/4; ie. bi aligned with the minor axis. The upper panel 
shows the angle k between the vectors 62 and 61 . The lattice 
configuration and basis vectors at times (a), (&), (c), (d) are 
shown in Fig. 4. 
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FIG. 4: The lattice configurations and the basis vectors in the 
fundamental region at different time (a), (6), (c), (d) indicated 
in Fig. 3. The discrete change of basis vectors between (b) 
and (c) is due to the mapping back to the fundamental region. 



